chr1 po1 gene1 chr2 po2 gene2 name
1 2 63456333 WDPCP 10 37493749 ANKRD30A TCGA-BH-A18M-01A
2 18 14563374 PARD6G 21 14995400 POTED TCGA-BH-A1ET-01A
3 10 37521495 ANKRD30A 3 49282645 CCDC36 TCGA-BH-A1ET-01A
4 10 37521495 ANKRD30A 7 100177212 LRCH4 TCGA-BH-A1EU-01A
5 18 18539803 ROCK1 18 112551 PARD6G TCGA-BH-A1EW-11B
6 12 4618159 C12orf4 18 1514414 PARD6G TCGA-C8-A12Q-01AWorkshop MOD: Visualisation de Données Omiques avec R
Introduction
Ce workshop MOD a été effectué par des étudiants du master Statistiques et Sciences des Données :
- Rim Alhajal (rim.alhajal@etu.umontpellier.fr)
- Abdoulaye Diop (abdoulaye.diop03@etu.umontpellier.fr)
- Pauline Dusfour-Castan (pauline.dusfour-castan@etu.umontpellier.fr)
- Lilou Zulewski (lilou.zulewski01@etu.umontpellier.fr)
Données
Les données utilisées dans ce Workshop sont disponibles dans la librairie OmicCircos de R. Elles reposent sur les bases de données fournies par The Cancer Genome Atlas Program (TCGA) contenant des données de patients atteints d’un cancer du sein. Voici la liste des tableaux utilisés :
TCGA.BC.cnv.2k.60 : Les variables de ce tableau sont respectivement le chromosome, sa position, le nom du gène observé ainsi que 60 autres variables représentant 60 individus différents.
TCGA.BC.fus : Ce tableau contient les données de fusion de gènes pour un individu précis.
TCGA.BC.gene.exp.2k.60 : Sur 60 individus, oune métrique liée à l’expression génique et le chromosome, sa position et le nom du gène observé sont consignés.
TCGA.BC.sample60 : Ce tableau référence le type de cancer (Normal, LumA, LumB, Basal ou Her2) dont sont atteints 60 personnes.
TCGA.BC_Her2_cnv_exp : Dans ce tableau, les t-value et p-value de tests statistiques et les métriques FDR et Bonferroni sont référencées.
TCGA.PAM50_genefu_hg18 : Ce tableau contient un ensemble de 50 gènes identifiés comme une signature d’expression génique associée aux sous-types du cancer du sein (Normal, LumA, LumB, Basal ou Her2) en fonction des schémas d’expression génique.
chr po Gene Basal Her2 LumA LumB Normal
1 1 212873846 CENPF 0.4829763 -0.02926616 -0.5437402 0.27822856 -0.07058307
2 1 120072157 PHGDH 0.6345189 -0.18662586 -0.3986822 -1.03013932 0.66043775
3 1 43599336 CDC20 0.3996952 0.00835552 -0.4690440 -0.07041247 -0.04550481
4 1 44992051 KIF2C 0.2035726 -0.16510205 -0.5053947 -0.18289071 -0.39001448
5 1 161575262 NUF2 0.4724002 -0.02381921 -0.7125208 0.58962688 -0.37053337
6 1 200572562 UBE2T 0.3899089 0.28453681 -0.5392594 0.73895213 -0.95238101Recueil des données
Afin de mettre en oeuvre les différentes méthodes de visualisation, il nous faut d’abord recueillir les données et les mettre en forme. Pour cela, suivez les étapes présentées ci-dessous :
- ÉTAPE 1 : Chargement de la librairie
Omiccircos
Pour pouvoir utiliser une librairie, il faut tout d’abord l’avoir installer à partir du CRAN puis l’appeler à l’aide de la fonction library().
library(OmicCircos)- ÉTAPE 2 : Importation des données
Nos données provenant du package OmicCircos (maintenant installé), nous les importons à l’aide de la fonction data().
data("TCGA.PAM50_genefu_hg18")
data("TCGA.BC.fus")
data("TCGA.BC.cnv.2k.60")
data("TCGA.BC.gene.exp.2k.60")
data("TCGA.BC.sample60")
data("TCGA.BC_Her2_cnv_exp")- ÉTAPE 3 : Mise en forme des données
# vecteur comprenant les positions de lignes des personnes
# possedant un cancer de type Her2
Her2.i = which(TCGA.BC.sample60[,2] == "Her2")
# vecteur comprenant les identifiants des personnes
# présentant un cancer de type Her2
Her2.n = TCGA.BC.sample60[Her2.i,1]
# indices des colonnes de TCGA.BC.cnv.2k.60 dont
# les identifiants sont présents dans Her2.n
Her2.j = which(colnames(TCGA.BC.cnv.2k.60)%in%Her2.n)
# creation d'un data d'expression des différents gènes
# pour chaque individus presentant un cancer Her2
cnv = TCGA.BC.cnv.2k.60[,c(1:3,Her2.j)]
# data cnv auquel on a enlevé les 3 premiere colonnes
cnv.m = cnv[,c(4:ncol(cnv))]
# bornage des valeurs >2 et <-2
cnv.m[cnv.m > 2] = 2
cnv.m[cnv.m < -2] = -2
# on ajoute aux 3 première colonne de cnv nos données bornées
cnv = cbind(cnv[,1:3], cnv.m)
# extraction des colonnes 1 à 3 de TCGA.BC.gene.exp.2k.60
# ainsi que des colonnes spécifiées par Her2.j
gene.exp = TCGA.BC.gene.exp.2k.60[,c(1:3,Her2.j)]Visualisations classiques
L’objectif de ces premières visualisations est d’assimiler la construction de graphes en reprenant les concepts de base de la librairie ggplot2 de R. Nous vous apprendrons notamment à tracer des graphes de densités et des histogrammes.
L’ensemble des visualisations de cette partie s’appuient sur la base de données TCGA.PAM50_genefu_hg18 et nécessite le chargement de divers packages.
# chargement des packages nécessaires
library(ggplot2)
library(dplyr)
library(tidyr)
library(plotly)
library(viridis)
library(viridisLite)
library(hrbrthemes)Les visualisations permettent d’étudier les différents types de cancer présentés dans le lexique#.
Visualisation des 10 premières observations des données
Voir le diagramme en barre
data(TCGA.PAM50_genefu_hg18)
df<-TCGA.PAM50_genefu_hg18
# Pour voir mieux, j'ai affiché que les 10 premières lignes du dataset. Si vous voulez voir plus, vous pouvez mettre en commentaire la ligne suivante.
df<-df[1:10,]
df_long <- df %>%
gather(key = "Subtype", value = "Expression", -chr, -po, -Gene)
ggplot(df_long, aes(x = Gene, y = Expression, fill = Subtype)) +
geom_bar(stat = "identity") +
labs(title = "Expression génique des gènes dans les sous-catégories de cancer du sein",
x = "Gène", y = "Expression génique") +
theme_minimal()
Utilisation de plotly
Plotly est une bibliothèque graphique interactive qui permet de créer des graphiques interactifs dans le langage de programmation R. Elle offre la possibilité de créer des graphiques tels que des diagrammes en barres, des diagrammes circulaires, des graphiques en nuage de points et bien d’autres, tout en permettant aux utilisateurs d’explorer et d’interagir avec les graphiques directement.
- Faisons un heatmap pour mieux voir les données dans l’ensemble
Voir le heatmap pour chaque type de cancer
df<-TCGA.PAM50_genefu_hg18
df_long <- df %>%
gather(key = "Subtype", value = "Expression", -chr, -po, -Gene)
# Créer le graphique ggplot
heatmap_plot <- ggplot(df_long, aes(x = Gene, y = Subtype, fill = Expression)) +
geom_tile() +
labs(title = "Heatmap de l'expression génique dans les sous-catégories de cancer du sein",
x = "Gène", y = "Sous-type de cancer du sein") +
theme_minimal()
# Convertir le graphique ggplot en heatmap interactif avec plotly
heatmap_interactif <- ggplotly(heatmap_plot, tooltip = "all")
# Afficher le heatmap interactif
heatmap_interactifUtilisation de ggplot2
Le package ggplot2 fonctionne par couches successives. La première d’entre elles, est un peu le canevas du graphe : elle consiste à indiquer dans quel jeu de données se trouve les données, et quelles sont les variables à représenter. Ensuite, une seconde couche est ajoutée, elle consiste à indiquer le type de graphe à réaliser : scatterplot, boxplot, barplot, etc… Enfin, les couches d’affinage permettent de choisir les couleurs, les échelles des axes, les options de légende et toutes autres options visuelles.
ggplotIl ne faut pas oublier le + entre les couches successives.
Comparaison de deux densités
Commençons par étudier certaines densités relatives aux types de cancer du sein deux à deux afin d’en retirer leurs différences.
Voici un exemple très basique de visualisation graphique de la comparaison des densités des types de cancer LumA et LumB qui sont ceux les plus communs :
# création de la base de données
data <- data.frame(value=c(TCGA.PAM50_genefu_hg18$LumA,
TCGA.PAM50_genefu_hg18$LumB),
type=c(rep("LumA", 50),
rep("LumB", 50)))
# création du graphe
comparaison <- data %>%
# création de l'objet
ggplot(aes(x=value, fill=type)) +
# définition du type de graphe
geom_density()
# affichage du graphe
print(comparaison)
%>%
Cet opérateur est un pipe, fréquemment représenté par une barre verticale | dans les terminaux : il renvoie la sortie d’une commande vers l’entrée d’une autre.
De nombreuses options de ggplot peuvent améliorer l’apparence générale des graphes. Dans l’exemple ci-dessous, vous pouvez retrouver respectivement la définission manuelle des couleurs, l’option de non-affichage de la légende associée à la coloration en fonction du type, ainsi que la personnalisation des éléments de texte et des légendes. Enfin, l’option alpha de geom_density règle la transparence de remplissage sous les courbes de densité.
# création du graphe
comparaison <- data %>%
# création de l'objet
ggplot(aes(x=value, fill=type)) +
# définition du type de graphe
geom_density(aes(color=type), alpha=0.4) +
# réglage des couleurs de légende
scale_color_manual(values=c("darkblue", "pink")) +
scale_fill_manual(values=c("darkblue", "pink")) +
# choix de l'affichage de la légende
guides(color="none") +
# personnalisation des éléments de texte
theme(strip.text.x=element_text(size=8),
plot.title=element_text(size=10, hjust=0.5, face="bold")) +
# titres
xlab("ARNm expression") +
ylab("Density") +
labs(title="Comparison of gene expression density across different types of cancer",
fill="Cancer type")
# affichage du graphe
print(comparaison)
À l’aide de ce graphe, nous pouvons remarquer que l’ARNmessager des gènes des deux types de cancer LumA et LumB ne s’expriment pas de la même façon. En effet, celui de LumA présente un pic autour de -0.5 tandis que celui de LumB présente un pic légèrement plus haut autour de 0.4.
Il est également possible de représenter les densités sur deux graphes différents en utilisant la commande facet_wrap(~XXX).
Utilisation de
facet_wrap()
Les échelles des différents graphiques peuvent être formatées différemment à l’aide de l’option scales : "fixed" permet d’obtenir des échelles égales pour tous les graphes mais d’autres options telles que "free_x", "free_y" ou "free" sont envisageables.
# création du graphe
comparaison <- data %>%
# création de l'objet
ggplot(aes(x=value, fill=type)) +
# définition du type de graphe
geom_density(aes(color=type), alpha=0.4) +
# séparation du graphe en sous parties
facet_wrap(~type, scales="fixed") +
# réglage des couleurs de légende
scale_color_manual(values=c("darkblue", "pink")) +
scale_fill_manual(values=c("darkblue", "pink")) +
# choix de l'affichage de la légende
guides(color="none") +
# personnalisation des éléments de texte
theme(strip.text.x=element_text(size = 8),
plot.title=element_text(size=10, hjust = 0.5, face="bold")) +
# titres
xlab("ARNm expression") +
ylab("Density") +
labs(title="Comparison of gene expression density across different types of cancer",
fill="Cancer type")
# affichage du graphe
print(comparaison)
À vous de jouer
En vous appuyant sur l’exemple précédent, remplissez les trous du code ci-dessous afin de représenter graphiquement la comparaison des densités densités d’expression d’ARNmessager pour les gènes des types de cancer Her2 et Normal.
Les données relatives aux cancer de type Normal sont contenues dans la variable Normal de la base de données.
Par la suite, les couleurs utilisées pour la représentation de Her2 et Normal seront respectivement darkred et darksalmon.
# création de la base de données
data <- data.frame(value=c(TCGA.PAM50_genefu_hg18$XXX,
TCGA.PAM50_genefu_hg18$XXX),
type=c(rep("XXX", 50),
rep("XXX", 50)))
# création du graphe
comparaison <- data %>%
# création de l'objet
ggplot(aes(x=XXX, fill=XXX)) +
# définition du type de graphe
geom_density(aes(color=XXX), alpha=0.4) +
# réglage des couleurs de légende
scale_color_manual(values=c("XXX", "XXX")) +
scale_fill_manual(values=c("XXX", "XXX")) +
# choix de l'affichage de la légende
guides(color="none") +
# personnalisation des éléments de texte
theme(strip.text.x=element_text(size=8),
plot.title=element_text(size=10, hjust=0.5, face="bold")) +
# titres
xlab("ARNm expression") +
ylab("Density") +
labs(title="Comparison of gene expression density across different types of cancer",
fill="Cancer type")
# affichage du graphe
print(comparaison)Afficher la Solution
# création de la base de données
data <- data.frame(value=c(TCGA.PAM50_genefu_hg18$Her2,
TCGA.PAM50_genefu_hg18$Normal),
type=c(rep("Her2", 50),
rep("Normal", 50)))
# création du graphe
comparaison <- data %>%
# création de l'objet
ggplot(aes(x=value, fill=type)) +
# définition du type de graphe
geom_density(aes(color=type), alpha=0.4) +
# réglage des couleurs de légende
scale_color_manual(values=c("darkred", "darksalmon")) +
scale_fill_manual(values=c("darkred", "darksalmon")) +
# choix de l'affichage de la légende
guides(color="none") +
# personnalisation des éléments de texte
theme(strip.text.x=element_text(size=8),
plot.title=element_text(size=10, hjust=0.5, face="bold")) +
# titres
xlab("ARNm expression") +
ylab("Density") +
labs(title="Comparison of gene expression density across different types of cancer",
fill="Cancer type")
# affichage du graphe
print(comparaison)
La densité de l’expression de l’ARNmessager de ce type de cancer est beaucoup plus étalée que celle du type de cancer Her2 et se mélange ainsi facilement dans les profils génomiques des autres types de cancer. Ce graphique met alors en évidence la rareté du type de cancer Normal et la difficulté à le classer.
Comparaison de toutes les densités
À vous de jouer
Pour mieux visualiser les différences entre tous les types de cancer, représentez la comparaison de la densité de l’ensemble des types de cancer.
Il est important de préciser au graphe que la variable doit être considéré comme une densité à l’aide de la commande aes(y=after_stat(density)).
Afficher la Solution
# création de la base de données
data <- data.frame(type=c(rep("LumA", 50),
rep("LumB", 50),
rep("Basal", 50),
rep("Her2", 50),
rep("Normal", 50)),
subtype=rep(TCGA.PAM50_genefu_hg18$chr, each=50),
value=c(TCGA.PAM50_genefu_hg18$LumA,
TCGA.PAM50_genefu_hg18$LumB,
TCGA.PAM50_genefu_hg18$Basal,
TCGA.PAM50_genefu_hg18$Her2,
TCGA.PAM50_genefu_hg18$Normal))
# densités des différents types de cancer
density <- data %>%
# création de l'objet
ggplot(aes(x=value, fill=type)) +
# définition du type de graphe
geom_density(aes(y=after_stat(density),color=type),
linewidth=0.5, alpha=0.4, fill=NA) +
# réglage des couleurs de légende
scale_color_manual(values=c("#69b3a2",
"darkred",
"darkblue",
"pink",
"darksalmon")) +
scale_fill_manual(values=c("#69b3a2",
"darkred",
"darkblue",
"pink",
"darksalmon")) +
# personnalisation des éléments de texte
theme(strip.text.x=element_text(size=8),
plot.title=element_text(size=10, hjust=0.5, face="bold")) +
# titres
xlab("ARNm expression") +
ylab("Density") +
labs(title="Comparison of gene expression density across different types of cancer",
color="Cancer type")
# affichage du graphe
print(density)
Nous remarquons à nouveau une densité d’expression de l’ARNm très étalée pour le type de cancer du sein Normal. De plus, les densités pour les types de cancer Basal et LumB sont assez similaires : Basal présente un pic légèrement plus haut et et une courbe plus concentrée autour de 0.5. Globalement, chacun des types de cancer possède une expression d’ARNm concentrée à différentes valeurs pour les gènes étudiés.
POUR ALLER PLUS LOIN…
Ajoutez les histogrammes respectifs des différents types de cancer sur le graphe précédent.
Afficher la Solution
# histogrammes et densités des différents types de cancer
density <- data %>%
# créationd de l'objet
ggplot(aes(x=value, fill =type)) +
# définition du premier type de graphe
geom_histogram(aes(y=after_stat(density)),
color="#e9ecef",
alpha=0.3,
binwidth=0.5,
position="dodge") +
# définition du second type de graphe
geom_density(aes(y=after_stat(density), color=type), linewidth=0.5, fill=NA) +
# réglage des couleurs de légende
scale_color_manual(values=c("#69b3a2", "darkred", "darkblue", "pink", "darksalmon")) +
scale_fill_manual(values=c("#69b3a2", "darkred", "darkblue", "pink", "darksalmon")) +
# choix de l'affichage de la légende
guides(color="none") +
# personnalisation des éléments de texte
theme(strip.text.x=element_text(size=8),
plot.title=element_text(size=10, hjust=0.5, face="bold")) +
# titres
xlab("ARNm expression") +
ylab("Density") +
labs(title="Comparison of gene expression density across different types of cancer",
fill="Cancer type")
# affichage du graphe
print(density)
geom_histogram
Cette commande du package ggplot2 permet de créer des histogrammes sur R. De nombreuses options sont disponibles comme color, alpha (règle la transparence du remplissage), binwidth (règle la largeur des intervalles) et position = "identity", "dodge" (règle la position des bandes).
Comparaison des densités pour chacun des chromosomes
À vous de jouer
En appliquant les compétences acquises jusqu’ici, remplissez les trous afin de représenter la densité et l’histogramme d’expression de l’ARNmessager des gènes pour le type de cancer Her2 pour chacun des chromosomes.
# création de la base de données
data <- data.frame(type=XXX,
value=XXX)
# densité pour chaque chromosome
density <- data %>%
# création de l'objet
ggplot(aes(x=XXX, bins=30, color=XXX, fill=XXX)) +
# définition du premier type de graphe
geom_histogram(alpha=0.6, position="identity") +
# définition du second type de graphe
geom_density(aes(y=after_stat(ndensity)),
alpha=0.5,
fill="grey",
color="black",
linewidth=0.25) +
# séparation du graphe en sous parties
facet_wrap(~XXX, scales="XXX", nrow=3) +
# réglage des couleurs de légende
scale_fill_viridis(discrete=FALSE) +
scale_color_viridis(discrete=FALSE) +
# choix de l'affichage de la légende
guides(color="none", fill="none") +
# personnalisation des éléments de texte
theme(strip.text.x=element_text(size=8),
plot.title=element_text(size=10, hjust=0.5, face="bold")) +
# titres
xlab("Chromosome") +
ylab("Frequency") +
labs(title="Gene mRNA expression densities for Her2 type of cancer on each of the chromosomes")
# affichage du graphe
print(density)Afficher la Solution
# création de la base de données
data <- data.frame(type=TCGA.PAM50_genefu_hg18$chr,
value=TCGA.PAM50_genefu_hg18$Her2)
# densité pour chaque chromosome
density <- data %>%
# création de l'objet
ggplot(aes(x=value, bins=30, color=type, fill=type)) +
# définition du premier type de graphe
geom_histogram(alpha=0.6, position="identity") +
# définition du second type de graphe
geom_density(aes(y=after_stat(ndensity)),
alpha=0.5,
fill="grey",
color="black",
linewidth=0.25) +
# séparation du graphe en sous parties
facet_wrap(~type, scales="fixed", nrow=3) +
# réglage des couleurs de légende
scale_fill_viridis(discrete=FALSE) +
scale_color_viridis(discrete=FALSE) +
# choix de l'affichage de la légende
guides(color="none", fill="none") +
# personnalisation des éléments de texte
theme(strip.text.x=element_text(size=8),
plot.title=element_text(size=10, hjust=0.5, face="bold")) +
# titres
xlab("Chromosome") +
ylab("Frequency") +
labs(title="Gene mRNA expression densities for Her2 type of cancer on each of the chromosomes")
# affichage du graphe
print(density)
Notons que la densité ne peut pas être représentée avec un seul point ce qui explique l’absence de courbe pour les chromosomes 3, 14 et 22. Quant à lui, l’histogramme représente le nombre de gènes ayant une valeur de test à un endroit donné : par exemple, pour le chromosome 16, il y a deux gènes qui ont une valeur de test environ égale à 0.21 pour le type de cancer Her2.
Pour la plupart des chromosomes, la densité observée est bimodale. Quelques chromosomes possèdent une densité différente. En effet, le chromosome 1 possède un grand pic sur l’histogramme autour de la valeur de test 0. Le chromosome 17 est aussi atypique puisque la densité associée est assez étalée dû aux 4 pics éloignés présents de l’histogramme.
POUR ALLER PLUS LOIN…
De la même façon que le graphe précédent, remplissez les trous pour représenter graphiquement la densité et les histogrammes d’expression de l’ARNm des gènes pour chaque type de cancer sur chacun des chromosomes.
# création de la base de données pour le graphe
data <- data.frame(type=TCGA.PAM50_genefu_hg18$XXX,
subtype=c(rep("LumA", 50),
rep("LumB", 50),
rep("Basal", 50),
rep("Her2", 50),
rep("Normal", 50)),
value=c(XXX,
XXX,
XXX,
XXX,
XXX))
# densité pour chaque chromosome
density <- data %>%
# création de l'objet
ggplot(aes(x=XXX, fill=XXX, color=XXX)) +
# définition du premier type de graphe
geom_histogram(alpha=0.6, position="dodge") +
# définition du second type de graphe
geom_density(aes(y=after_stat(ndensity), color=XXX),
alpha=0.25,
fill="grey",
linewidth=0.25) +
# séparation du graphe en sous parties
facet_wrap(~XXX, scales="XXX", nrow=3) +
# réglage des couleurs de légende
scale_color_manual(values=c("#69b3a2",
"darkred",
"darkblue",
"pink",
"darksalmon")) +
scale_fill_manual(values=c("#69b3a2",
"darkred",
"darkblue",
"pink",
"darksalmon")) +
# choix de l'affichage de la légende
guides(color="none") +
# personnalisation des éléments de texte
theme(strip.text.x=element_text(size=8),
plot.title=element_text(size=10, hjust=0.5, face="bold")) +
# titres
xlab("Chromosome") +
ylab("Frequency") +
labs(title="Gene mRNA expression densities for each type of cancer on each of the chromosomes",
fill="Cancer type")
# affichage du graphe
print(density)Afficher la Solution
# création de la base de données pour le graphe
data <- data.frame(type=TCGA.PAM50_genefu_hg18$chr,
subtype=c(rep("LumA", 50),
rep("LumB", 50),
rep("Basal", 50),
rep("Her2", 50),
rep("Normal", 50)),
value=c(TCGA.PAM50_genefu_hg18$LumA,
TCGA.PAM50_genefu_hg18$LumB,
TCGA.PAM50_genefu_hg18$Basal,
TCGA.PAM50_genefu_hg18$Her2,
TCGA.PAM50_genefu_hg18$Normal))
# densité pour chaque chromosome
density <- data %>%
# création de l'objet
ggplot(aes(x=value, fill=subtype, color=subtype)) +
# définition du premier type de graphe
geom_histogram(alpha=0.6, position="dodge") +
# définition du second type de graphe
geom_density(aes(y=after_stat(ndensity), color=subtype),
alpha=0.25,
fill="grey",
linewidth=0.25) +
# séparation du graphe en sous parties
facet_wrap(~type, scales="fixed", nrow=3) +
# réglage des couleurs de légende
scale_color_manual(values=c("#69b3a2",
"darkred",
"darkblue",
"pink",
"darksalmon")) +
scale_fill_manual(values=c("#69b3a2",
"darkred",
"darkblue",
"pink",
"darksalmon")) +
# choix de l'affichage de la légende
guides(color="none") +
# personnalisation des éléments de texte
theme(strip.text.x=element_text(size=8),
plot.title=element_text(size=10, hjust=0.5, face="bold")) +
# titres
xlab("Chromosome") +
ylab("Frequency") +
labs(title="Gene mRNA expression densities for each type of cancer on each of the chromosomes",
fill="Cancer type")
# affichage du graphe
print(density)
En comparant l’expression de l’ARNm pour l’ensemble des types de cancer, nous remarquons que les conclusions faites pour le type de cancer Her2 peuvent être utilisées à nouveau. En effet, les densités associées sont généralement bimodales sauf piur quelques chromosomes. De plus, nous identifions une expression d’ARNm plus concentrée sur les chromosomes 8, 16 et 17 pour le type de cancer Her2 ainsi que sur les chromosomes 7, 11 et 12 pour le type de cancer LumA. Ce grahe permet de mettre en lumière les différentes expressions d’ARNm de mêmes gènes pour les différents types de cancer.
Création de HEATMAP
L’objectif de cette partie est de réaliser des heatmap intéractives représentant les observations.
La librairie heatmaply dont nous allons nous servir est incluse dans la librairie plotly (plotly) qui est une librairie permettant d’inclure de l’interactivité dans les graphiques.
Dans notre exemple de heatmap, vous pourrez notament connaître la valeur d’une cellule en faisant simplement glisser votre souris dessus, zommer sur sur les cellules et se déplacer dans la figure.
Pour cette première visualisation, il vous faut donc charger la librairie “heatmaply” et faire un choix de gradient de couleurs à l’aide des commandes suivantes :
library(heatmaply)
gradient_col <- ggplot2::scale_fill_gradient2(
low = "blue", high = "red", midpoint = 0, limits = c(-8, 8))Vous pouvez trouver ci-joint deux sites donnant des exemples d’utilisation et de l’aide sur la librairie heatmaply :
Dans un premier temps, nous allons sélectionner le chromosome que l’on souhaite étudier.
Prenons ici le chromosome numéro 1.
# on selectionne le chromosome souhaité dans notre dataframe
chr1 <- subset(gene.exp, gene.exp$chr==1)
# on selectionne les colonnes qui nous interessent
mydata_1 <- chr1[, 3:18]
# nom des lignes
rownames(mydata_1) <- mydata_1[, 1]
# on supprime la première colonne qui ne nous interesse plus
mydata_1 <- mydata_1[, -1] Création de la heatmap à l’aide de heatmaply :
heatmaply(mydata_1, scale_fill_gradient_fun = gradient_col,
Colv=NA, Rowv=NA, scale='none', limits = c(-8, 8),
xlab = "Individus",ylab = "Gènes",
main = "Expression des gènes du chromosome 1")La heatmap ci-dessus nous donne visuellement une idée de l’expression de chaque gène du chromosome 1. On observe que dans la majorité des cas le gène “CLCA2” est sur-exprimé (les cellules apparaissent très rouges), et les gènes “PDZK1” “DLANI1” sont sous-exprimés (les cellules apparaissent plutôt bleutées).
Vous pouvez choisir d’afficher les corrélations avec la fonction heatmaply_cor() ou les valeurs manquantes avec la fonction heatmaply_na().
A vous de jouer
En vous basant sur l’exemple ci-dessus, créez une heatmap représentant l’expression des gènes du chromosome 17.
Voici une idée du résultat que vous devez obtenir :
Afficher la Solution
chr17 <- subset(gene.exp, gene.exp$chr == 17 )
mydata_17 <- chr17[,3:18]
rownames(mydata_17) <- mydata_17[,1]
mydata_17 <- mydata_17[,-1]
heatmaply(mydata_17, scale_fill_gradient_fun = gradient_col,
Colv=NA, Rowv=NA, scale='none', limits = c(-7, 8),
xlab = "Individus",ylab = "Gènes",
main = "Expression des gènes du chromosome 17")Création d’une heatmap circulaire
L’objectif de cette seconde partie sur les heatmap est de réaliser des heatmap circulaires avec dendrogrammes.
Commencez par charger les packages nécessaires :
library(circlize)
library(ComplexHeatmap)On effectue une petite modification de notre dataframe en supprimant une colonne qui ne nous sera pas utile par la suite.
# suppression de la colonne position
gene.exp2 <- gene.exp[, -2] Pour pouvoir réaliser notre heatmap, il nous faut tout d’abord convertir notre dataframe en matrice :
mat <- as.matrix(gene.exp2, labels=TRUE) # on transforme le dataframe en matrice
# On définit les noms des colonnes et des lignes
row_names <- as.vector(gene.exp2$NAME)
col_names <- c("chr", "NAME", "TCGA.A2.A04W.01A", "TCGA.A2.A0D1.01A",
"TCGA.A2.A0EQ.01A", "TCGA.A8.A081.01A", "TCGA.A8.A08J.01A",
"TCGA.A8.A09X.01A", "TCGA.AN.A0FV.01A", "TCGA.AO.A0J2.01A",
"TCGA.AO.A12D.01A", "TCGA.B6.A0RS.01A", "TCGA.BH.A0EE.01A",
"TCGA.C8.A12L.01A", "TCGA.C8.A12P.01A", "TCGA.D8.A13Z.01A",
"TCGA.E2.A1B0.01A")
rownames_factor <- factor(row_names)
# On met notre matrice sous format "numeric"
mat_numeric <- matrix(as.numeric(mat), ncol = ncol(mat))
rownames(mat_numeric) = row_names
colnames(mat_numeric) = col_names
# Vecteur des chromosomes :
split <- as.vector(gene.exp2$chr)
split <- factor(split)Le vecteur splitdoit contenir les différents groupes/catégories si vous souhaitez divisez votre heatmap par groupes/catégories.
Création des heatmap circulaires :
- Sur cette première heatmap, nous avons affiché tous les chromosomes ainsi que tous les gènes de la manière la plus basique possible.
## Affichage des heatmaps
circos.clear()
circos.par(start.degree = 90, gap.degree = 3)
col_fun1 = colorRamp2(c(-6, 0, 6), c("blue", "white", "red"))
# heatmap simple
circos.heatmap(mat_numeric[, !(colnames(mat_numeric) %in% c("chr", "NAME"))],
split = split,
col = col_fun1)
- Sur cette seconde heatmap, nous avons encadré chaque chromosome et affiché son numéro.
# heatmap avec les chromosomes associes
circos.clear()
circos.heatmap(mat_numeric[, !(colnames(mat_numeric) %in% c("chr", "NAME"))],
split = split, col = col_fun1, track.height = 0.4,
bg.border = "green", bg.lwd = 2, bg.lty = 2,
show.sector.labels = TRUE)
- Sur cette troisième heatmap, nous avons rajouté le nom de chaque gène.
# heatmap avec les chromosomes et les noms des genes
circos.clear()
circos.heatmap(mat_numeric[, !(colnames(mat_numeric) %in% c("chr", "NAME"))],
split = split, col = col_fun1, track.height = 0.4,
bg.border = "green", bg.lwd = 2, bg.lty = 2,
show.sector.labels = TRUE,
rownames.side = "outside", rownames.cex = 0.4)
- Sur cette quatrième heatmap, nous avons rajouté les dendrogrammes ainsi qu’une légende de couleur.
Note
La méthode de clustering appliquée est hclust soit une méthode de classiification hierarchique, la distance utilisée par défaut est celle de Ward. Vous pouvez trouver des informations sur hclustau Lien suivant et sur la classification hierarchique ici et ici.
# heatmap avec les dendrogrammes
circos.clear()
circos.heatmap(mat_numeric[, !(colnames(mat_numeric) %in% c("chr", "NAME"))],
split = split, col = col_fun1, track.height = 0.4,
bg.border = "green", bg.lwd = 2, bg.lty = 2,
show.sector.labels = TRUE,
dend.side = "outside")
lgd = Legend(title = "values", col_fun = col_fun1)
grid.draw(lgd)
Le nombre de chromosomes étant important, il est compliqué de voir en détail les gènes. L’exercice suivant consiste donc à sélectionner certains chromosomes et à afficher la heatmap circulaire correspondante.
Il est important de toujours appeler la commande circos.clear() avant chaque nouvelle heatmap.
A vous de jouer
Sélectionnez les chromosomes 1, 3, 6, 11, 17 et 19.
x <- c("XXX") # vecteur des chromosomes que l'on souhaite selectionner
data_chr <- gene.exp2[gene.exp2$chr %in% x, ] # dataframe avec nos chromosomes selectionnes
# Convertion en matrice
mat_chr <- as.matrix(data_chr, labels=TRUE)
mat_chr_numeric <- matrix(as.numeric(mat_chr), ncol = ncol(mat_chr))
# On définit les noms des colonnes et des lignes
row_names_chr <- as.vector(data_chr$NAME)
col_names_chr <- c("chr", "NAME", "TCGA.A2.A04W.01A", "TCGA.A2.A0D1.01A",
"TCGA.A2.A0EQ.01A", "TCGA.A8.A081.01A", "TCGA.A8.A08J.01A",
"TCGA.A8.A09X.01A", "TCGA.AN.A0FV.01A", "TCGA.AO.A0J2.01A",
"TCGA.AO.A12D.01A", "TCGA.B6.A0RS.01A", "TCGA.BH.A0EE.01A",
"TCGA.C8.A12L.01A", "TCGA.C8.A12P.01A", "TCGA.D8.A13Z.01A",
"TCGA.E2.A1B0.01A")
rownames(mat_chr_numeric) = row_names_chr
colnames(mat_chr_numeric) = col_names_chr
# Vecteur des chromosomes :
split2 <- as.vector(data_chr$chr)
split2 <- factor(split2)Réaliser maintenant une heatmap affichant les noms des gènes ainsi que les dendrogrammes.
Voici une idée du résultat que vous devez obtenir :
Afficher la Solution
x <- c(1, 3, 6, 11, 17, 19)
data_chr <- gene.exp2[gene.exp2$chr %in% x, ]
mat_chr <- as.matrix(data_chr, labels=TRUE)
mat_chr_numeric <- matrix(as.numeric(mat_chr), ncol = ncol(mat_chr))
row_names_chr <- as.vector(data_chr$NAME)
col_names_chr <- c("chr", "NAME", "TCGA.A2.A04W.01A", "TCGA.A2.A0D1.01A",
"TCGA.A2.A0EQ.01A", "TCGA.A8.A081.01A", "TCGA.A8.A08J.01A",
"TCGA.A8.A09X.01A", "TCGA.AN.A0FV.01A", "TCGA.AO.A0J2.01A",
"TCGA.AO.A12D.01A", "TCGA.B6.A0RS.01A", "TCGA.BH.A0EE.01A",
"TCGA.C8.A12L.01A", "TCGA.C8.A12P.01A", "TCGA.D8.A13Z.01A",
"TCGA.E2.A1B0.01A")
rownames(mat_chr_numeric) = row_names_chr
colnames(mat_chr_numeric) = col_names_chr
split2 <- as.vector(data_chr$chr)
split2 <- factor(split2)
circos.clear()
circos.heatmap(mat_chr_numeric[, !(colnames(mat_chr_numeric) %in% c("chr", "NAME"))],
split = split2, col = col_fun1, track.height = 0.4,
bg.border = "green", bg.lwd = 2, bg.lty = 2,
show.sector.labels = TRUE, dend.side = "outside",
rownames.side = "inside", rownames.cex = 0.35)
lgd = Legend(title = "values", col_fun = col_fun1)
grid.draw(lgd)
Vous pouvez consulter le lien suivant si vous souhaitez en apprendre plus sur les heatmap circulaires : heatmap circulaire
Visualisation avec Circos
A présent, nous allons réaliser une visualisation circulaire à l’aide du package omicCircos et de la fonctioncircos.
Vous pouvez trouver de l’aide au package OmicCircos aux liens suivants :
colors = rainbow(10, alpha = 0.5) # choix des couleurs
# Si vous souhaitez enregistrer votre image sous forme de pdf
#pdf("visucomp_histo.pdf", 8,8) # nom de l'image
#par(mar = c(2, 2, 2, 2)) # definition de la fenetre graphique
plot(c(1,800), c(1,800), type = "n", axes = FALSE, xlab = "", ylab = "", main = "")
zoom = c(1, 22, 939245.5, 154143883, 0, 180)
###### Moitié droite du cercle ##########
# mettre les chromosomes et l'emplacement des genes
circos(type = "chr", R = 400, cir = "hg18", W = 4,
print.chr.lab = TRUE, scale = TRUE, zoom = zoom)
# creation de la heatmap
circos(type = "heatmap2", R = 300, cir = "hg18", W = 100,
mapping = gene.exp, col.v = 4, cluster = FALSE,
col.bar = TRUE, lwd = 0.01, zoom = zoom)
# nombre de copies (gain_rouge, perte_bleu)
circos(type = "ml3", R = 220, cir = "hg18", W = 80, mapping = cnv,
col.v = 4, B = FALSE, lwd = 1, cutoff = 0,
zoom = zoom)
# creation de l'histogramme
circos(type = "h", R = 140, cir = "hg18", W = 80, mapping = cnv,
col.v = 4, B = TRUE, lwd = 1, col = colors[1],
zoom = zoom)
# lien entre certains genes
circos(type = "link", R = 130, cir = "hg18", W = 10,
mapping = TCGA.BC.fus, lwd = 2, zoom = zoom)
# Mise en valeur des chromosomes 11 et 17
the.col1 = rainbow(10, alpha = 0.5)[1]
highlight = c(140, 400, 11, 282412.5, 11, 133770314.5, the.col1, the.col1)
circos(type = "hl", R = 110, cir = "hg18", W = 40,
mapping = highlight, lwd = 2, zoom = zoom)
the.col2 = rainbow(10, alpha = 0.5)[6];
highlight = c(140, 400, 17, 739525, 17, 78385909, the.col2, the.col2)
circos(type = "hl", R = 110, cir = "hg18", W = 40,
mapping = highlight, lwd = 2, zoom = zoom)
highlight.link1 = c(400, 400, 140, 376.8544, 384.0021,
450, 540.5)
circos(type = "highlight.link", cir = "hg18",
mapping = highlight.link1, col = the.col1, lwd = 1)
highlight.link2 = c(400, 400, 140, 419.1154, 423.3032,
543, 627)
circos(type = "highlight.link", cir = "hg18",
mapping = highlight.link2, col = the.col2, lwd = 1)
####### Moitié gauche du cercle #######
# Chromosome 11
zoom = c(11, 11, 282412.5, 133770314.5, 180, 270)
circos(type = "chr",R = 400, cir = "hg18", W = 4, print.chr.lab = TRUE,
scale = TRUE, zoom = zoom)
circos(type = "heatmap2", R = 300, cir = "hg18", W = 100,
mapping = gene.exp, col.v = 4, cluster = FALSE, lwd = 0.01,
zoom = zoom)
circos(type = "ml3", R = 220, cir = "hg18", W = 80, mapping = cnv,
col.v = 4, B = FALSE, lwd = 1, cutoff = 0,
zoom = zoom)
circos(type = "h", R = 140, cir = "hg18", W = 80, mapping = cnv,
col.v = 4, B = TRUE, lwd = 1, col = colors[1],
zoom = zoom)
# Chromosome 17
zoom = c(17, 17, 739525, 78385909, 274, 356)
circos(type = "chr", R = 400, cir = "hg18", W = 4, print.chr.lab = TRUE,
scale = TRUE, zoom = zoom)
circos(type = "heatmap2", R = 300, cir = "hg18", W = 100,
mapping = gene.exp, col.v = 4, cluster = FALSE, lwd = 0.01,
zoom = zoom)
circos(type = "ml3", R = 220, cir = "hg18", W = 80, mapping = cnv,
col.v = 4, B = FALSE, lwd = 1, cutoff = 0,
zoom = zoom)
circos(type = "h", R = 140, cir = "hg18", W = 80, mapping = cnv,
col.v = 4, B = TRUE, lwd = 1, col = colors[1],
zoom = zoom)
De l’extérieur vers l’intérieur, vous pouvez trouver :
le numéro du chromosome et la position des gènes
une heatmap de l’expression des gènes
le nombre de copies de chaque gènes par rapport à la normale (rouge : gain, bleu : perte)
des histogrammes pour observer la dispersion de la loi de distribution
des liaisons entre les protéines de fusion (tracées avec l’algorithme de courbe de Bézier)
Lexique
hg18 : La représentation du génome humain a connu plusieurs versions avant d’aboutir à la représentation la plus optimale. La version 18 correspond à l’avant dernière version.
Her2 : Il s’agit d’une protéine naturellement présente dans l’organisme : un récepteur transmembranaire impliqué dans la régulation de la prolifération cellulaire.
LumA : Le cancer du sein Luminal A est l’un des sous-types luminaux et est généralement associé à des caractéristiques moins agressives que Luminal B.
LumB : Ce type de cancer du sein est associé à un grade plus élevé, des taux de prolifération accrus, et un pronostic global plus défavorable.
Basal : Il s’agit d’un sous-type du cancer du sein négatif pour les récepteurs d’oestrogène (ER), les récepteurs de progestérone (PR) et le récepteur 2 du facteur de croissance épidermique humain (HER2). On le désigne souvent comme un cancer du sein triple négatif (TNBC).
FDR (First Division Restitution) : Ceci réfère à l’évènement où l’une des cellules filles produites lors de la première division méiotique conserve les deux chromatides d’un chromosome, sans subir la séparation normale en deux cellules distinctes. Cela aboutit à un gamète avec un ensemble complet de chromosomes, plutôt que la moitié attendue.
Bonferroni : Cette correction statistique permet d’ajuster le niveau de confiance de chaque intervalle pour que le niveau de confiance simultané ne soit pas erroné.